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O . Abstract 

O 



The retarded Green's function G(k, lu) of a single small spin polaron in the 



j^ ■ three-band model for the Cu02 plane is calculated in the self-consistent Born 

approximation. It is shown that such a spin polaron is a good quasiparticle 
excitation for realistic values of spin exchange J and effective hopping r. The 
polaron spectral density Ap{]i.,uj) demonstrates small damping in contrast to 
the results of calculations starting from the bare hole, i.e. the pole strength of 
the energetically low-lying quasiparticle peak ^p(k) varies from 50% to 82% 



for J JT oc 0.1 -^ 0.7. The quasiparticle peak dispersion reproduces the main 
features of the bare polaron spectrum 17^ near the band bottom. 

The sphericaUy symmetric approach is used for the description of spin 
excitations. It makes it possible to consider the quantum antiferromagnetic 
background without the spontaneous symmetry breaking and the unit ceh 
doubhng. 

The new method of the self-consistent calculation, based on continued 

fraction expansion of Green's function, is represented in details. The method 

preserves the proper analytical properties of the Green's function and provides 

the possibility to analyze the nature of its singularities. 
71.27.+a,71.10.Fd,75.30.Mb 
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I. INTRODUCTION 

Up to date much theoretical work has been devoted to the problem of a hole motion in two 
dimensional (2D) s=l/2 quantum antiferromagnet (AFM) |jl[]. The important question is 
whether a hole injected in the undoped ground state behaves like a quasiparticle (QP). This 
problem is mainly investigated within the framework of self-consistent Born approximation 
(SCBA) for the t — J model [@-0] and Kondo-lattice [^. There are rather few works devoted 
to the three-band Hubbard model or the Emery model [^,Tn| that is more realistic for Cu02 



planes in high- T^ superconductors (HTSC). For the t — J model it was shown that the 
spectral density function Ah(k, u) of a doped hole revealed a QP peak of intensity Zk ~ J/t 
and a broad incoherent part that has a width of about 6 — 7t. The QP band bottom 
corresponds to the momenta ki = (±7r/2, ±7r/2) . Similar results were obtained for the 



Emery model ||TT|JT^ . The presence of a great incoherent part and small intensity of QP 



peak indicate, that bare holes are rather poor elementary excitations even for k close to ki 

In order to investigate the hole motion in the t — J model one usually decouple the hole 
operator into a spinless fermion and an antiferromagnetic magnon operator. As a result the 
zero approximation corresponds to the dispersionless band with zero energy of the hole. The 
hopping of the particle appears only due to the fermion-magnon scattering which is treated 
by the usual perturbation method in k-space. For this reason we think that in this approach 
the resulting quasiparticle pole in the fermion Green's function takes into account mainly a 
polaron of a large radius. The close situation takes place in the usual treatment of a hole 
motion in the effective three-band model [p!l| -|T3| and the Kondo-lattice model [§], when one 
starts from a bare hole, not from a magnetic polaron of a small radius. 

In the present paper, within the framework of the effective three-band model we study the 
spectral density function Ap(k, u) of a single small polaron, i.e., an excitation which takes into 
account from the very beginning a local hole-spin coupling . It is well known that the simplest 
candidate for such a small polaron is an analog of the so-called Zhang- Rice singlet in CUO4 



plaquette |T^, |]T5[. The mean- field spectrum fik of tliis excitation is well studied |T^ and 



it will be used as the zero approximation in our treatment. We shall consider the coupling 
of small polaron to spin-wave excitations within the self-consistent Born approximation 
(SCBA) for the corresponding two-time retarded Green's function G{\i.,uj). 

Our motivations to study Ap(k,u!) and the corresponding quasiparticle band are the 
following. First, it may be easily shown for the one- hole problem that the mean- field energy 
of the polaron i7k represents the center of gravity of the spectral function 

/oo 
ujAp{k,uj)diu (1) 

-oo 

This means that the minimum of fl]^ is the upper bound of the actual position for the 
quasiparticle band bottom. The SCBA, based on a bare hole Green's function, gives the 



minimum value of quasiparticle energy equal to ujh,mm = — 2.6r [|12|, for a typical value 
of copper-copper AFM exchange constant J = 0.7r. Here r is a constant of the effective 
oxygen-oxygen hopping via an intervening copper site (note that our unit of energy is twice 
that of Ref. |T2[ r = 2t). As to the value of the small polaron mean field band bottom it 
turns out to be substantially lower than ujh,min, ^k = — 3.17r, for the same value of J/r. So 
we may conclude that important local correlations are lost in SCBA when we start from the 
bare hole operators. 

Secondly, we shall show that a small polaron represents the elementary hole excitation 
much more better than a bare hole dressed by magnons within the framework of SCBA. This 
will be manifested by a relatively great intensity of a quasiparticle peak in our calculation. 

Finally, the mean field spectrum fi^ of the simplest small spin polaron explicitly depends 
on the state of the antiferromagnetic background . In the case of long range order state, 
r2k demonstrates a flat band region close to the magnetic Brillouin zone boundary [jl5 . 



This region corresponds to the bottom of the band . Moreover, if one takes into account 
the direct oxygen-oxygen hopping, flnite temperature and a more complicated form for a 
small polaron wave function, then fik reproduces the experimentally observed extended 
saddle point [|T6|-pO[] which is directed along the line (0, vr) — (0,0) [^. Therefore it seems 



important to clear up if the quasiparticle band reproduces the pecuharities of fl^ dispersion. 
Below we study this question for the simplest variant of the model. 

The distinctive feature of our investigation consists in considering the AFM copper spin 
subsystem in a spherically symmetric approach [^,^]. Such an approach is the most ap- 



propriate for treating the quantum 2D AFM at any finite temperature. As a result the 
scattering of a spin polaron by spin excitations in the singlet spin background leads to the 
spectral function periodicity relative to the full Brillouin zone. Note that the conventional 
two-sublattice spin approach leads to periodicity relative to the magnetic (reduced) Brillouin 
zone ||-0p|. 



The paper is organized as follows. In Section II we give the derivations for the self- 
consistent equation for the Green's function in the case of the small polaron approach. In 
Section III we present the procedure that provides the opportunity to avoid the iterative 
solution of the self-consistent equation for complex energies. The procedure is based on the 
continued fraction (CF) expansion of Green's function and it gives a possibility to calculate 
consequently the coefficients of the CF expansion with the use of the quadrature method. 
In Section IV we offer the termination of the continued fraction which leads to the correct 
analytical properties of the resulting Green's function. Numerical results for the self-energies 
and spectral functions, the relation of our results to the previous approaches and discussion 
are given in Section V. Section VI summarizes the results. An Appendix contains some 
details of the approach which gives the expression for integrals over the spectral density in 
terms of the chain representation of continuous fraction. 

Part of our results were presented in Brief Reports p4| . The paper presented here gives 



further results as well as the description of the new method and more details about the 
calculations. 



II. EFFECTIVE HAMILTONIAN AND SMALL POLARON GREEN'S 

FUNCTION 

Following Refs. [^,TD,T3| we adopt the Hainiltonian that corresponds to one hole problem 
in the Cu02 plane in HTSC: 

r,ai,a2,(T,(T' i",g 

ai, a2 = ±lgx, ±lgy, g = ±gx, ±gy 

Here and below g^^y are basic vectors of a copper square lattice ( |g| = 1), r + a are four 
vectors of O sites nearest to the Cu site r, the operator cj. creates a hole with the spin index 
cr = ±1 at the O site, Sa-a' = \(}aa' , operator S represents the localized spin on the copper 
site. As mentioned above r is the integral of oxygen hole hoppings that takes into account 
the coupling of the hole motion with copper spin subsystem. J is the constant of nearest 
neighbor AFM exchange between the copper spins. 

It is well known that the most prominent feature of the Hamiltonian (j^) is that the 
low-energy physics of hole excitations is described by the Bloch sums B]^ ^ based on one site 
small polaron operators Sj^ 

eL-^Ee-BL. (3) 

^ a 
r,r' 

Here {,},[,] stand for an anticommutator and commutator, respectively; (...) = 
Q^^Sple''^^ ...], Q = Spe~^^ . Sp implies taking the trace of an operator, j3 = {kT)~^ is 
an inverse temperature; a = —a; Z^^'^^ = |o"i) (o"2| are the Hubbard projection operators for 
Cu sites states; 7k = :| X]gexp(«k ■ g). 

To calculate the average for commutators and anticommutators such as K]^, we take 
into account that these expressions are reduced to the two-site or three-site averages of the 

6 



Hubbard operators located at different copper sites. In the spherically symmetric approach 
adopted here, all these averages can be expressed in terms of the two-site spin correlation 
functions Cr = (So ■ Sr). Let us recall that due to the spherical symmetry (SfSj) = 
Sap (S^Sj") = |(Sj-Sj), (S'f) = 0. Simultaneously, at T ^ the spin subsystem is 
described by a long range order state with finite effective magnetization m, Cr(|r| -^ oo) = 
rn^l^—lY^+^y ^ the value of m is dictated by the Bose condensation of spin excitations at the 
antiferromagnetic vector qo = (7r,7r). 

Note, that Bl^^\L) corresponds to the Cu02 plane state with the total spin equal to | if 
\L) is the singlet state. We treat i3]J ^ as a candidate for the elementary excitations operator 
and calculate the corresponding two-time retarded Green's function (^(k, uj) and spectral 
density 

Ap{k, uj) = ^^(k, UJ + iO+), 

TT 
/•oo 

Gik,u) = {B^,Ma)^ = -^ / rfte^"(*-*')({^k,.(t),4.(0})- (5) 

Jt' 

The Dyson's equation for G(k, u) has the form 

G-\k,uj) = G,' - S(k,.;),S(k,c.) = (7^k,.|<,)(^'-^^ (6) 



where 



Go = (cu - f]k)-^7^k,. = [B^,^,H] = ^==^e-^''■'•7^r,., (7) 



7^r,<. = -4rSr,. + K, + ni,, (8) 



\g,a,o-i g, a, 0-1,0-2 

^ \g,a,oi 



Cr+g+a,o-2 I ! l"J 



0201 N 

r ;'^r+a,oi I ) 



fik = ({7^k,o, < J) = {rQ, + JQj)lK^. (10) 

Q.(k) = -^ - 8 (^ + ^g) 7k + (^ - Q + \C^^ 72k + 2 (i - Cg + iCd) 7dk, 



Qj(k)=Cg(7k-4), 

Here and below d = g^, + gy, jdk = cos{k^g) cos{kyg). 

We see from Eqs.(^) and (pH]), that the self-energy E(k, u;) accounting for interaction ef- 
fects is expressed through the higher-order Green's functions. One should notice, first, that 
the terms linear in i3k,cr do not contribute to the irreducible Green's function ([Tl|). Second, 
the lowest-order self-energy contribution is provided by the first term in the right-hand side of 
the expression (^TJ), while the second term leads to higher-order corrections. Following Ref. 
we evaluate (|ll]) with a proper decoupling procedure for two time correlation function 
('^k,o-(i^)'^ko-(^'))- This procedure is equivalent to the self-consistent Born approximation 
in a usual diagrammatic technique 0. In our case this means that the two-time correla- 
tion function is decoupled into the spin-spin correlation function and the polaron-polaron 
correlation function. The adopted decoupling procedure preserves the main character of 
polaron site operator (^) - four hole site operators surround the copper spin operator. It 
schematically has the form 



Let us mention that we also tested the more complex decoupling procedure, and it did not 
qualitatively alter the results given by approximation (p!2D . On the next step we project 
polaron operators in (|l^) onto B^a '■ 

c,{t)Z^{t)c^^BUt),^=({c^{t)Z,it),Bl^}). (13) 



Taking into account, that we shall calculate only the irreducible part of Green's function 
(|n|), the averages {Zy.j^(t)Zy.^{t')) are transformed to corresponding spin-spin correlation 
functions (5*° (t)S'"^(t')V Collecting all terms, we have 

8 



/Clc- 



(^k,.(tX,(t'))^Ar-i^^r2(k,q)(Sk-,,.(t)4-,,.(t'))(S-c(t)S,(t')), (14) 



K, 



J, 



r(k,q) = rr.(k,q) + -rj(k,q), 
r,(k, q) = 47k-q [(1 + 7k-q) /2irk-q - 1] 

(s_,(t)s,(t')) = ^E^^"-^^"""^ (Sr(t)sxt')) 

r,r' 

Using the spectral representation for Green's functions we obtain the following interme- 
diate result for the self-energy 



E(k,.;)=iV-iE^r2(k,q)/ 



°° duji 



■X 



oo TT 



dw2 ePi^^+^^) + 1 Im[G'(k - q, c^i -t- zO+)]Im[D(q, cj2 + «0+)] 
-oo TT (e^'^i + 1) (e^'^2 _ 1) ^ _ (^^ + ^2) + «0+ 

Here the spin excitation Green's function p2| , p3| 

8JCg l-7q . 



(15) 



D(q,u;) = (SlJ5^) 



3 C<j2 — CJq ' 



(16) 



OOr, 



-32JaiCg(l-7q)(2A + l + 7q) 



we neglect the influence of doped holes on copper spin dynamics and take the spin spectrum 
parameters calculated in Ref. [^] (the vertex correction ai = 1.7 , the spin excitations 
condensation part m? = 0.0225). 

As a result we come to the integral equation for the Green's function that always arises 
within the framework of SCBA 



G(k,cu) 



1 



00 — ^k — S(k, cu) 



(17) 



where 



S(k, uj) = N-' E M\k, q) [(1 + z/q)G(k -ci,u-uj^) + u^G{k - q, cu + u^)] , (18) 



z/q = 1/ [exp(/3c<jq) — 1] is the Bose function, 



r2/U „N _ ^k-q^2„, „^ (-4Cg) (1 - 7q; 



M^(k,q) = ^r^(k,q) ^ -^^'^ '^' . (19) 

r(k, q) corresponds to the bare vertex for the couphng between a spin polaron and a 
spin wave. It is known [^ that this vertex is substantially renormalized for q close to the 
AFM vector qg = (tt, vr) . This renormalization is due to the strong interaction of a polaron 
with the condensation part of spin excitations that must be taken into account from the 



very beginning. As a result, the renormalized vertex r(k, q) must be proportional to ||25[ 



[q — qo) + L^^ , Lg being the spin-spin correlation length, L^ ^ oo in our case of a long 
range order state of the spin subsystem. Below this renormalization is taken into account 
empirically by the substitution 



r(k,q) ^ r(k,q) = r(k,q)^(l + 7q) . (20) 

The introduced vertex correction is proportional to |q — qo| for q close to qo . We have used 



also the following two functions for the vertex correction Jl + 7^ and w 1 + 7^ and have 
obtained the results similar to those presented below. Let us mention that the bare vertex 
leads the dramatic decrease of the QP bandwidth. 

III. SOLUTION OF THE INTEGRAL EQUATION 



The equation ([T7| ) is usually solved by an iteration procedure. We propose here an 
alternative way, based on the continued fraction expansion of (^(k, z) 

G(k,^) = — ^° ^i ^^ , a„ = a„(k), 6„ = 6„(k), (21) 

Z — CLq — z — di — z — d^ — 

where 

/+00 ][ /■+00 

Ap{\s.,Lo)duj = K]^, '^0 = To ujAp(k,Lu)duj = fik- 

The coefficients 6„, o,„, n > are related with the spectral density Ap(k, uo) via the set of 
orthogonal polynomials P„(u;), satisfying the recurrence 

P_i(cu)=0, Pq{uj) = 1, 
10 



Pn+l{uj) = {UJ - an)Pn{uj) - blPn^l{uj), 



(22) 



and 



bl^i 



J^::P^{u;)A,{k,u;)du;' 



J^^P^^,{u;)A,{k,uj)duj 



(23) 



(24) 



Jl^P^iuj)A,ik,uj)du; 
Here we have used the nonnormahzed form of the polynomials 

J^^PnHPsHA,ik,uj)duj = dnsimZi bmf- 



Comparing Eqs.(^) and (|T^ we see that the self energy S(k, z) is the continued fraction 
similar to ^(k, z) . Thus we can introduce the spectral density 



p{k,uj) = -Im[E(k,cu + zO+)]/7r 



and the set of polynomials n„(ci;) with the recurrence analogous to (^21): 



where 



0"n+l 



n„(u;) = (cL^-a„)n„_i(u;) -6^n„_2(t^), no{uj) = 1, n_i(t^) = 0, 

/oo ]^ r+oo 

p(k,uj)du, '^1 = 72 / ci;p(k, ci;)^^^, 
-oo Oi J —oo 






n>l. 



On the other hand, we have from Eq.(18 



(25) 



p(k, uj) = — iV-i ^ M2(k, q) [(1 + z/q)Ap(k - q, ^ - a;q) + v^Ap{k - q, ^ + u^)] . (26) 

-f^k q 



If we shall put the expression for p(k, cu) in Eq.(p5D then the coefficients a^+i and 6„+2 
will be expressed through the integrals of the form 



/oo 
(u;±u;q)™n,^y4p(k-q,u;)rfw, i < n, m = 0, 1, 
-oo 



(27) 



Now, the trick is that the polynomials in uj in the integrals ( P7D have the degree less 
than or equal 2n+ 1. As it was proved by Nex ^\ such integrals may be expressed through 
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the coefficients {ao, . . . , a„, 60 • • • , &n}-The details of such a procedure are presented in the 
Appendix. So, it turns out that in SCBA we can recursively calculate pairs of coefficients 
On+i, 6„+iand to obtain S(k, uj) in the continued fraction form. Of course we must calculate 
the coefficients at all the chosen k + q points in the ffist Brillouin zone simultaneously. 
Below the chosen points correspond to a lattice of 32 x 32 unit cells. The presented procedure 
provides the opportunity to avoid the iterative solution of Eq. ([TtD for complex energies. 

IV. TERMINATION OF THE CONTINUED FRACTION 

The procedure outlined in the previous section would be efficient if after calculating a 
finite number of coefficients an,bn,n < uq, we could appropriately approximate the part 
(infinite in our case) of the continued fraction T„p , which has not been calculated. In other 
words, we rewrite the expression ( pTD in the form 

G(k, .) = -^ '-^ \ ,, . , (28) 

z-ao-z-ai- z - ttno - Tn^,{k, z) 

and try to find a function T„p (so called "terminator") that is close to T„(,. 

Various ways to construct such approximations are described in the literature on the 
recursion method (see Refs. |p7H29(] ). The asymptotic behavior of continued fraction co- 



efficients is governed by the band structure and singularities of spectral density p9|. The 
main asymptotic behavior depends on the band structure: {a„} and {6„} converge towards 
limits in the single band case, oscillate endlessly in a predictable way in the multiband case. 
Damped oscillations are created by isolated singularities. The main point here is that an 
isolated simple pole produces exponentially daraped contribution in {a^}, {bn} , n — > 00. For 
our case it means that the quasiparticle pole position and weight could be obtained with the 
high accuracy from finite number of coefficients, and the asymptotic behavior determines 
the incoherent part of the spectrum. It is obvious that the spectrum we deal with has an 
lower bound and no upper bound. We can thus expect that coefficients will not converge to 
some finite values but will grow up to infinity. 
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In Fig.|l] we represent the coefficients a„, bn as functions on n calculated according the 
procedure described in the preceding section for two values oi J, J = 0.7r and J = O.lr 
and for k = (7r/2,7r/2). As it may be seen, the distinctive feature of this dependence is 
that for large n the coefficients a„, 6„ are linear functions of n. For all this, the slope for 
a„ coefficients is twice as many as the slope for 6„. So the coefficients behavior may be 
approximated as 

bn ~ Xin + A2, a„ ^ 2\in + A3, A, = A,(k), n > 1 (29) 

It is interesting that the analogous behavior show the coefficients for t — J model approx- 
imated by slave-fermion Hamiltonian p[ treated within the self-consistent Born approxi- 
mation ^. For J = 0.4t, k = {11/2,71/2) the coefficients a„,6„ dictated by the relation 
analogous to (|18D are represented in Fig. (a). 

Now we shall show that the same asymptotic ( J29|) has the CF expansion of incomplete 
Gamma function that is written as |^ 

r(a,x) = ^:^l^^^^... (30) 

^ ' ^ a;+ 1+ x+ 1+ ^ ^ 

This circumstance we shall use for the construction of the terminator T/v(k, 2) for G(k, z) 



Let us introduce the function 

T(a, —x) ^ , ( 1 — a 



g\a,x) 



Q^{—xY 

where 



i/(^-r^)- (^'' 



In order to rewrite the CF ( pT]) in the form analogous to Eq.(^) we denote 

1 



Then it is not difficult to obtain the relations 
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1 + rit„. (33) 



" x-(2n + l-a) -(n + l)(n + l-a)t„+i' ^ ^ 

so that 

g{a,x)=to (35) 



has the form ( ]28|) with the coefficients % = 1 and 

On = 2n + 1 — a, 6^ = ?2(?2 — a). (36) 



Comparing Eqs. ( pQ]) and (|36|) for large n, when Jn{n — a) ^ n — a/2, we see that the 

substitution 

2A2 2; + 2A2 - A3 + Ai 

" = -v" = a; 

leads to the function G 

n(\ \ ^ ~i ^-^2 2 + 2A2 — A3 + Ai 

<^(k, ^) = );-^(-— ' 3;^ ) (37) 



that has the same asymptotic as G(k, 2;) (^ll) . It means that G'(k, 2;) can be used as the 
terminator for G(k, z). , i.e. we can express T„Q(k, z) = &„o+it„o+i through G'(k, 2;), and 
the coefficients a„, 6„, n < no, and then substitute it for r„o(k, z) (see [^ for the details of 
matching Greenians). 

We thus obtain G(k, z) in the whole complex energy plane including the real axis. Note 
that usually the procedure of discretizing the energy range lo is used for iteration process 
when the Dyson's equation is solved numerically. It is not obvious that such a self-consistent 
solution leads to the correct analytical properties of the resulting Green's function. In 
contrast the continued fraction representation guarantees these properties (e.g. the positive 
definiteness of spectral function Ap(k, u;)). 

V. RESULTS AND DISCUSSION 

In this section we represent our results for the retarded Green's function (^(k, tu) for 
three-band model at T = 0. The self consistent equation (|I^) was solved on a 32 x 32 cell 
lattice. The number of calculated CF levels no was taken no = 30. 
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To begin with, we check the vahdity of the method outhned above by calculating the 
spinless hole Green's function for the t — J model and compare the results with the results 
of Martinez and Horsch [Q obtained by usual iteration procedure. In Fig.^ y4/j(ki,cij + 
^r]), ReS(ki,u;), — ImS(ki,co'), ki = (7r/2,7r/2) for the value of J = 0.4t are represented 
for the 16x16 site lattice and broadening constant t] = O.Olt. The comparison of Fig.|^ 
(b),(c),(d) and the corresponding functions given by Fig. 7,8 from Ref. 0] (the same lattice 
size and t]) demonstrates that the position of peaks of the hole spectral function and their 
intensities coincide. Some difference is that our Ah{\<.i,uj), is more smooth and there are no 
strong oscillations in our self energy S(k, cu) in the interval (— 2t < uj < — 0.75t). 

The results for the small spin polaron spectral density, real and imaginary parts of the 
self energy for the characteristic value of energetic parameter J = 0.7r are given in Fig. |^ 
for the symmetrical points ki = (7r/2, 7r/2), k2 = (0, 0), ks = (vr, vr). The energy broadening 
parameter is ?7 = 0.002 (we will refer all quantities in units of r from now on). We find main 
common feature in the spectral density for ki and k2: a sharp quasiparticle peak exists at 
the bottom of each spectrum. The position of the QP peak corresponds to the condition 
ReG~^(k, cu) = , i.e., the point where we have the crossing of the functions y = u — Q]^ 
and ReS(k, cu) , see Figs.^ (a) and ^ (b). In Fig. ^ (d) we show Ap(ki,u) calculated for 
rj = 0.002 (solid line) and rj = 0.0005 (dashed line) in order to study the scaling behavior of 
the peaks and their widths with respect to changes in rj. Both peaks fit quite closely with 
a Lorentzian (l/vr) lZ{]ii)T]/ (to — e(ki)) + r]'^ > , e(ki) is the location of the peak, which 
in the limit rj —^ becomes Z(\<.i)6{uj — e(ki)). It means that ImS(k, cjp) — i> in the same 
limit. Here and below we shall speak about the position e(k) of such peaks (with imaginary 
part of pole close to zero) as about the energy of the quasiparticle. 

Figs. (a) and ^ (b) also demonstrate that the incoherent part of Ap(k, u) increases and 
the pole strength decreases with the increase of e(k), Z(ki) = 0.82, Z(k2) = 0.347. Let us 
recall that i7k represents the center of gravity of the spectral function. In our Figures the 
center of gravity corresponds to the crossing of the real axis by the line y = uj — Q]f_. So, if 
the quasiparticle peak is far from this point we would have a great incoherent part. 
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Quite different features demonstrate Ap(k3,u) in Fig.^ (c). The broad lowest peak is 
determined by the appearance of nonzero ImS(k3,u;) in the region where ReG~^(k3,u;) has 
no zeroes. Two broad additional peaks at uj ^ —1.6 and uj ~ —1.05 are formed due to 
the zero values of Re(j>'~^(k3, tu) close to these u. But at the same time the ImZl(k3,u;) is 
strong in these regions. Moreover, the maximum of ImE(k3,a;) (near the point u ^ —1.37 
) determines the local minimum of Ap{]<.s,uj) despite of the fact, that this point is close to 
the frequency value where ReG'~^(k3,u;) = 0. It is clear that it is impossible to treat any 
of Ap{]<.s,uj) peaks as a quasiparticle one. Let us mind that the qualitative behavior of the 
real part of the self-energy in Fig.|^ (c) is close to that one which is represented by Kampf 
and Schrieffer, see ||3^ Fig.^] (b), for the pseudogap regime of the Hubbard model. Fig.|^ 



(c) demonstrate three solutions of ReG~^(k3,C(j) = 0. Although there is a sharp crossover 
from a situation with three solutions to one quasiparticle solution, the spectral function still 
changes smoothly due to the presence of imaginary part of E. 

Figs. (b) and (c) demonstrate qualitatively different character of Ap(k, u) for the points 
k2 = (0,0), k3 = (7r,7r). This is the consequence of our spherically symmetric approach for 
treating the AFM copper spin subsystem. As mentioned in Introduction this approach gives 
rise to the spectral function periodicity relative to the full Brillouin zone, not magnetic one. 

In Fig.§ we show the dispersion relation e(k) of the quasiparticle band and mean field 
dispersion fi^ along symmetry lines in the Brillouin zone. For e(k) we reproduce only that 
k— values for which the lowest peak has a pronounced quasiparticle peak taking the following 
criteria: -ImS(k, e(k.) + tri) < 21], rj = 0.002. As it is known [|15|, due to the antiferromagnetic 
character of spin correlation functions the Q\^ demonstrates a "fiat dispersion region" close 
to the line 7k < 0, |7k| <^ 1 ,i.e. close to the boundary of magnetic Brillouin zone X-N-X, 
see Fig. ^. As it may be seen from Fig. |^ the quasiparticle band exists for greater part of 
the Brillouin zone except the region of the top of the f^k spectrum. Moreover, the dispersion 
law e(k) qualitatively reproduces main features of the spectrum fik. As it is mentioned in 
the Introduction, fik demonstrates the important features of the hole spectrum for Cu02 



plane if one takes into account 0-0 hoppings and spin frustration |£T|. We hope that in this 
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case e(k) will reproduce these features also. 

Let us compare the small polaron spectral density ylp(k, cj) with the results for the bare 
hole Ah{k,uj) given by Kabanov and Vagov |T^ for ki = (7r/2,7r/2), J = 0.7r , see Fig. |^ 
(a). First, Fig. ^ (a) indicates that Ap(ki, u) demonstrates much sharper QP peaks relative 
to the results for a bare hole. For example, the pole strength Zp(ki) for the QP peak of 
y4p(k, u) equals to Zp{k.) = 0.82 . The corresponding value for Ah given by Ref. |jl2| is much 
smaller, ^^(k) = 0.25. 

Secondly, Fig.|| (a) explicitly demonstrates the one peak structure Ap(ki, u) in contrast to 
Ah{\i.i,uj). Finally, it is important, that the bottom of our QP e(k) = —3.52 is substantially 
lower than u^inA = —2.6 from [|12[. These results are the consequence of the fact, that 



elementary excitation - spin polaron of small radii iSk,^ - from the beginning takes into 
account the strong local hole-spin coupling. 

It is clear that the QP peaks for a bare hole and a small polaron have to coincide in the 
exact solution of the problem. The mentioned above discrepancies between our calculations 
and that one of Ref. ^^ are the consequence of different approximations. 



In order to test the convergence of our results relative to the increase of the lattice size 
and the number of calculated CF levels NL, in Fig.^ we show the QP peak of Ap(\<.,uj) at 
k = (7r/2, 7t/2) , J = 0.7 for different lattices and hq. This peak, as it is seen from Fig.|] (a) 
and (b), changes insignificantly in going from 24 x 24 to 32 x 32 cell lattice and from 22 to 
30 NL. 

We consider now the transformation of G{\<., u) with the decrease of J. In order to clarify 
how the character of Ap(k,uj) peaks is changed, in Fig.^ we show Ap(k,u!) for the value of 
J = 0.1 at points ki = (7r/2, 7r/2), k2 = (0,0), ks = {11,71). The decrease of J leads to the 
enlargement of broad incoherent part of y4p(k, cj). 

As before the fiat band region of the QP band bottom enlarges along a magnetic Brillouin 
zone boundary. It is represented by the point ki . In Fig.^ (a) ^^(ki,^^) demonstrates 
explicitly a rather strong QP peak, Zp(ki) = 0.5 , which corresponds to the condition 
ReG~^(k,u!) = 0. Quite different character of Ap(k,u) is typical for k that correspond to 
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the tops of flk band: in the low-energy sector for k2,k3, see Fig.H (b),(c), one observes 
Ap{k,uj) - peaks with small intensity. For example, the pole strength Zp of such a QP peak 
for Ap{'k2,uj) equals to Zp{k2) = 0.016. Taking uji{k.) as the value of uj corresponding to the 
center of these lowest by energy peaks one can see that ReG'~^(k, co'i(k) ^ for the k under 
discussion. Figs. ^ (b) and (c) demonstrates that these peaks are determined by the peaks 
in ImE(k, C(j) at points c<j;(k). The self-energy part S(k, c<j) occurs through Green's function 
of a small polaron bounded to a spin waves. These peaks can be considered as the QP band 
of such complex states. 

If we shall treat the width of the quasiparticle band W as the difference uJi{]i.2 = (0,0)) — 
a;i(ki = (7r/2,7r/2)) then W turns out to be of order J for small values of J (J ~ 0.1) in 
line with the results for the hole Green's function approach [|12|. 

It is clear, that for small J/t the concept of a small spin polaron fails and it is important 
to estimate the validity limits of this concept. Our calculations demonstrate, that the inten- 
sity of QP peaks and the structure v4p(k, uj) do not change dramatically for k, corresponding 
to the band bottom, up to J/t = 0.1. For example, ^p(f , |) ~ Zp(n, 0) ~ 0.50 at J/r = 0.1. 
So the J/t lowest boundary value of the small spin polaron concept validity is lower than 
J/t = 0.1. 

In Table 1 we give the numerical values u!i{\<.) of the center position of the lowest Ap{]<.2, uj) 
peaks (a;;(k) = e(k) for k values where QP peak is observed) and their pole strength (area 
under the peak) Zp(k) for k = {n/2, tt/2), (0, 0), (0, tt) and different values of J. 

We do not represent the results for large J, J ^ t, as our approach in the present form 
fails to describe this limit. Here, from the very beginning we treat a small polaron by a single 
site operator i?r ( ^. For large J the mean field static energy of such a state is proportional 
to J and such a state is unstable. So in this limit we must extend the basis of site operators. 
The simplest way is to include in the basis the additional operator of a bare hole. In SCBA 
approximation this will lead to the system of two selfconsistent equations. As a result all 
effects of spin subsystem - hole interaction will be proportional to t/J. The more general 
procedure for the extending of the small polaron operator basis is outlined in pl |. 



VI. SUMMARY 

We have studied the small spin polaron motion in the three-band model. The two- 
time retarded Green's function was calculated within the framework of self-consistent Born 
approximation for 32 x 32 cell lattice. We have shown that spin polaron of small radius 
represents a good approximation to the true quasiparticle low-energy excitation even at 
mean-field level. The account of the self-energy does not crucially change the polaron motion 
picture for realistic values of parameters. For quasimomenta k values, corresponding to the 
band bottom, most of the total spectral weight is concentrated in the quasiparticle peak 
(Table 1). In the same region of k-space the shape of the QP dispersion curve e(k) reproduces 
that of mean- field dispersion fik(Fig.^. 



We compare our results with previous investigations ||Tl| , pr2| that start from the bare 
hole. We see that the small polaron mean-field energy fik lies much lower than the QP pole 
obtained from SCBA for the bare hole. Since Q^ determines the center of gravity for the 
Green's function spectral density the actual QP pole position (at least for the band bottom) 
should lie deeper in energy than f2k (Figs. ^ and ^.This means that in the three-band 
model the important local correlations should be taken into account in zero approximation 
and small spin polaron should be constructed. Then the polaron scattering on spin waves 
is of less importance and it may be treated by perturbation methods. 

The conclusion is that the low-energy physics of high- Tc superconductors should be 
considered in terms of small spin polaron dynamics. In particular, that the problem of 
superconducting hole pairing must be treated as pairing of these quasiparticles rather than 
pairing of bare holes. 
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APPENDIX. CHAIN REPRESENTATION AND INTEGRALS OVER THE 

SPECTRAL DENSITY 



The integration over spectral density could be done with the quadrature approach |2^. It 



is very efficient one when applied to the electron structure calculations |^. Unfortunately, 
the spectral density we deal with has no upper bound and exponentially depends on the 
energy, i.e. it substantially differs from typical spectral density that appears in a band- 
structure calculations. It turned out that the direct application of Nex's quadrature approach 
p7| is not stable numerically for our purposes. That is why for calculations of the integrals 
(^) we use the chain representation of continuous fraction. It means that the CF of the form 



(^) may be interpreted as the Green's function G{uj) = (uo\{uj — h)^^\uo) of the one particle 
tight-binding Hamiltonian h of the semi-infinite one-dimensional lattice, the a„, bn, \un) being 
the site energies, nearest neighbor hoppings and on-site basis states respectively 



^n+l 



an = \Un\h\Un) , &n+l = (un\h\Ur, 

We introduce the eigenstates \ipm), and eigenenergies Em of the chain Hamiltonian 

m 

then the Green's function spectral density becomes the local density of states on the zeroth 



site of the chain |E6 



A{uj) = lmG{uj + lO^) = V {uoli'm) 5iuj - Em) ii'mluo) . 

71 

m 

Then the following identities hold 

/ + 00 /• + 00 

f{uj)A{uj)duj= f{uj)J2{^o\i^m)S{uj - Em) {^pm\uo)duj = (38) 

-CO J —OO yy^ 
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Yl (Uol-ipm) f{Em) (^m|Mo) = (uo\f{h)\Uo 



(39) 



Nex has proved [^] that for a polynomial / of the degree 2no + 1 the integral F for the 
infinite chain has the same value as the analogous integral for the truncated chain of the 
length riQ + 1. The Hamiltonian of the truncated chain in the basis of states {\uq) . . . \uno)} 
has the form of the tri diagonal {no + 1) x [uq + 1) matrix 

ao hi 
hi ai 62 



a„„_i 6, 



no 



"no '-^no 



Now, instead of integrating the spectral density function by uo we directly calculate the 
matrix function f{hT) in order to take /oo matrix element. Then F = foo , as it follows 



from the last identity of Eq. (pQl). We see that the answer for F is expressed only through 
the first coefficients {qq, . . . , a„, 60 • • • , &n}- 
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FIGURES 
FIG. 1. The coefficients o„ and 6,„ of the continued fraction expansion of Gp{k,uj) as functions 

on n for k = (7r/2,-7r/2): a) J = 0.7r, b) J = O.lr. Calculated on the 32 x 32 cell lattice. 

FIG. 2. Results for the hole Green's function Gh{k,uj) for the t — J model calculated with the 
same parameters as in ref.( [Q]), J = 0.4t, k = (7r/2,7r/2), r/ = 0.01, 16 x 16 site lattice, a) The 
coefficients o„ and 6„ of the continued fraction expansion of Gh{k, u) as functions on n. b) Spectral 
function Ah{k,uj). c) Real part of the self-energy, d) Imaginary part of the self-energy. The unit 
of energy is t = 1. 

FIG. 3. Spin polaron spectral density Ap{k,uj), real and imaginary parts of the self-energy 
T,{k,uj) calculated for J = 0.7r and 32 x 32 cell lattice and different values of k. a) k = (7r/2,7r/2), 
here we also reproduce the hole spectral function Afi{k,io) obtained in Ref. [l3] and pointed with 
h-arrows; h) k = (0,0); c) k = (vr,7r). In Figures (a-c) rj = 0.002r, the sloping straight lines repre- 
sent the function m — r2(A;). d) The dependence of the quasiparticle peak of Ap{k = (7r/2,7r/2),a;) 
for two values of the broadening factor r]. The unit of energy is r = 1. 

FIG. 4. The dispersion of the QP band €{k) and the mean field dispersion fifc along the 
symmetry lines in the Brillouin zone (see inset) for J = 0.7, 32 x 32 cell lattice, rj = 0.002. 

FIG. 5. Ap{k = {7r/2,TT/2),Lo) for J = 0.7 and rj = 0.002 calculated for a) NL = 30 and 
different lattice sizes; b) 32 x 32 cell lattice and different numbers of calculated CF levels NL. 

FIG. 6. Spin polaron spectral density Ap{k,uj), real and imaginary parts of the self-energy 
T,{k,uj) calculated for J = O.lr, rj = 0.002, 32x32 cell lattice at three values of A;: a) A: = (7r/2,7r/2), 
h) k = (0,0), c) k = (7r,7r). 
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TABLES 
TABLE L Position of the lowest by energy peak uji{\i.) and the area under the peak ^p(k) for 

different values of J/t and k 



J/t 


Zp(0,0) 


c^Ko,o) 


Zp(7r/2,7r 


/2) 


(jK7r/2,7r/2) 


Zp{7T,0) 


COl{7T,0) 


0.1 


0.016 


-4.24 


0.50 




-4.48 


0.55 


-4.51 


0.3 


0.039 


-3.37 


0.72 




-4.09 


0.714 


-4.13 


0.5 


0.174 


-2.714 


0.793 




-3.79 


0.738 


-3.83 


0.7 


0.347 


-2.25 


0.823 




-3.52 


0.808 


-3.56 
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